
# Packages ----------------------------------------------------------------

# using a reproducible environment
renv::restore()

# pacman allows to check/install/load packages with a single call
# if (!require("pacman")) install.packages("pacman") # already in renv.lock
library("pacman")

# packages to load (and install if needed) -------------------------------
pacman::p_load(
  here,      # easy file paths
  see,       # theme_modern and okabeito palette
  report,    # reporting various info 
  labelled,  # labelled data
  # ---- packages specific to this project ----
  scales,
  lme4,
  simr,
  readxl,
  openxlsx,
  statmod,
  qqplotr,
  emmeans,
  ggbeeswarm,
  ggpubr,
  patchwork,
  quarto,
  easystats,
  # ---
  tidyverse  # modern R ecosystem
)


# Custom functions shared across scripts ----------------------------------
source(here("scripts/_functions.R"))


# Global cosmetic theme ---------------------------------------------------

theme_set(theme_modern(base_size = 14)) # from see in easystats

# setting my favourite palettes as ggplot2 defaults
options( 
  ggplot2.discrete.colour   = scale_colour_okabeito,
  ggplot2.discrete.fill     = scale_fill_okabeito,
  ggplot2.continuous.colour = scale_colour_viridis_c,
  ggplot2.continuous.fill   = scale_fill_viridis_c
)


# Fixing a seed for reproducibility ---------------------------------------
set.seed(14051998)


# Adding all packages' citations to a .bib --------------------------------
knitr::write_bib(c(.packages()), file = here("bibliography/packages.bib"))

1 Rationale

To account for the non-normal, positively skewed distributions of the RTs, we fitted Generalized Linear Mixed Models (GLMMs) with inverse Gaussian distributions. The models were implemented in the lme4 R package and integrated in tidymodels workflows using the package multilevelmod. Models with Gamma and Gaussian distributions were also fitted and compared with the AIC and BIC to ensure that we chose the best distribution available.

The models included the Group (aphantasic, control), Congruence condition (congruent or incongruent) and Color condition (color or uncolored) along with all their two and three way interactions as fixed categorical predictors, while participants have been included as grouping factors (i.e. “random effects”). The random effect structure was chosen by fitting and comparing models with every possible combination of distribution and structure (intercept by participant, congruence or color, slope by participant on congruence and/or color) aiming for the best balance between goodness of fit and parsimony. Complex random-effects structures including various slopes on the factors failed to converge to stable and reliable estimates, hence the optimal models chosen included a single by-participant random intercept.

2 Model fitting

The formula of the model fitted is \(RT \sim (aphantasia + congruence + color)^3 + (1|subjectid)\).

model_formula <- rt ~ (aphantasia + congruence + color)^3 + (1|subjectid)

fit_implicit <- 
  glmer(
    formula = model_formula,
    data = df_i_rt,
    family = inverse.gaussian(link = "identity"),
    control = glmerControl(optimizer = "bobyqa")
  )

fit_explicit <- 
  glmer(
    formula = model_formula,
    data = df_e_rt,
    family = inverse.gaussian(link = "identity"),
    control = glmerControl(optimizer = "bobyqa")
  )

3 Model diagnostics

The quality checks of the models are displayed in Figure S3.1 and Figure S3.2 .

Code

# characteristics to check
model_checks <- c("pp_check","homogeneity", "vif", "outliers", "qq", "reqq")

check_model(
  fit_implicit, 
  checks = model_checks, 
  detrend = FALSE, 
  residual_type = "normal"
  ) |> 
  plot()

ggsave(here("figures/checks-implicit.png"), dpi = 600)

Figure S3.1: Model assumption checks for the Generalized Linear Mixed Model fit on the RTs in the implicit task.

Code

check_model(
  fit_explicit, 
  checks = model_checks, 
  detrend = FALSE, 
  residual_type = "normal"
  ) |> 
  plot()

ggsave(here("figures/checks-explicit.png"), dpi = 600)

Figure S3.2: Model assumption checks for the Generalized Linear Mixed Model fit on the RTs in the explicit task.

4 Model summaries

Table S4.1: Performance and estimates of the GLMM fitted on the implicit task data.

AIC	AICc	BIC	R2 (cond.)	R2 (marg.)	ICC	RMSE	Sigma
-4619.06	-4619.03	-4548.26	0.01	3.47e-04	0.01	0.26	0.48

# Fixed Effects
Parameter	Coefficient	SE	95% CI	t(8767)	p
(Intercept)	0.67	0.01	(0.64, 0.69)	53.99	< .001
aphantasia1	0.02	0.02	(-0.03, 0.07)	0.77	0.442
congruence1	-0.01	4.03e-03	(-0.02, -5.69e-03)	-3.37	< .001
color1	-2.57e-03	4.02e-03	(-0.01, 5.32e-03)	-0.64	0.523
aphantasia1 × congruence1	0.02	8.05e-03	(7.70e-03, 0.04)	2.92	0.004
aphantasia1 × color1	-7.96e-04	8.04e-03	(-0.02, 0.01)	-0.10	0.921
congruence1 × color1	0.01	8.03e-03	(-5.65e-03, 0.03)	1.26	0.209
aphantasia1 × congruence1 × color1	0.02	0.02	(-0.01, 0.05)	1.20	0.229

# Random Effects
Parameter	Coefficient
SD (Intercept: subjectid)	0.08
SD (Residual)	0.48

Table S4.2: Performance and estimates of the GLMM fitted on the explicit task data.

AIC	AICc	BIC	R2 (cond.)	R2 (marg.)	ICC	RMSE	Sigma
-2673.61	-2673.58	-2603.01	0.03	7.96e-04	0.03	0.28	0.41

# Fixed Effects
Parameter	Coefficient	SE	95% CI	t(8593)	p
(Intercept)	0.79	0.02	(0.75, 0.83)	37.14	< .001
aphantasia1	0.01	0.04	(-0.07, 0.09)	0.25	0.799
congruence1	-0.02	4.21e-03	(-0.03, -9.82e-03)	-4.29	< .001
color1	-0.03	4.20e-03	(-0.04, -0.02)	-6.93	< .001
aphantasia1 × congruence1	0.03	8.42e-03	(0.01, 0.04)	3.19	0.001
aphantasia1 × color1	7.19e-03	8.39e-03	(-9.26e-03, 0.02)	0.86	0.392
congruence1 × color1	2.74e-03	8.39e-03	(-0.01, 0.02)	0.33	0.744
aphantasia1 × congruence1 × color1	0.01	0.02	(-0.02, 0.05)	0.84	0.402

# Random Effects
Parameter	Coefficient
SD (Intercept: subjectid)	0.12
SD (Residual)	0.41

5 Estimated means and contrasts

Table S5.1: Estimated means for each group in each congruence condition and contrasts between groups and conditions in the implicit task.

Group	Condition	Median (ms)	SE	df	asymp.LCL	asymp.UCL
Control	Incongruent	669.71	19.36	Inf	631.77	707.65
Aphantasia	Incongruent	676.89	15.84	Inf	645.84	707.93
Control	Congruent	644.38	19.31	Inf	606.53	682.22
Aphantasia	Congruent	675.05	15.83	Inf	644.02	706.08

Marginal Contrasts Analysis
Level1	Level2	Difference	95% CI	SE	df	z	p
Aphantasia Incongruent	Aphantasia Congruent	1.84e-03	(-0.01, 0.01)	5.23e-03	Inf	0.35	0.725
Aphantasia Incongruent	Control Congruent	0.03	(-0.02, 0.08)	0.02	Inf	1.30	0.193
Control Congruent	Aphantasia Congruent	-0.03	(-0.08, 0.02)	0.02	Inf	-1.23	0.219
Control Incongruent	Aphantasia Congruent	-5.34e-03	(-0.05, 0.04)	0.02	Inf	-0.21	0.831
Control Incongruent	Aphantasia Incongruent	-7.18e-03	(-0.06, 0.04)	0.02	Inf	-0.29	0.774
Control Incongruent	Control Congruent	0.03	( 0.01, 0.04)	6.13e-03	Inf	4.13	< .001

Marginal contrasts estimated at aphantasia, congruence p-values are uncorrected.

Table S5.2: Estimated means for each group in each congruence condition and contrasts between groups and conditions in the explicit task.

Group	Condition	Median (ms)	SE	df	asymp.LCL	asymp.UCL
Control	Incongruent	801.78	33.22	Inf	736.68	866.88
Aphantasia	Incongruent	799.21	27.05	Inf	746.18	852.23
Control	Congruent	770.29	33.17	Inf	705.28	835.30
Aphantasia	Congruent	794.57	27.04	Inf	741.56	847.57

Marginal Contrasts Analysis
Level1	Level2	Difference	95% CI	SE	df	z	p
Aphantasia Incongruent	Aphantasia Congruent	4.64e-03	(-0.01, 0.02)	5.45e-03	Inf	0.85	0.394
Aphantasia Incongruent	Control Congruent	0.03	(-0.05, 0.11)	0.04	Inf	0.68	0.499
Control Congruent	Aphantasia Congruent	-0.02	(-0.11, 0.06)	0.04	Inf	-0.57	0.570
Control Incongruent	Aphantasia Congruent	7.21e-03	(-0.08, 0.09)	0.04	Inf	0.17	0.866
Control Incongruent	Aphantasia Incongruent	2.57e-03	(-0.08, 0.09)	0.04	Inf	0.06	0.952
Control Incongruent	Control Congruent	0.03	( 0.02, 0.04)	6.42e-03	Inf	4.91	< .001

Marginal contrasts estimated at aphantasia, congruence p-values are uncorrected.

6 Visualisations

The figures below are also displayed in the main article. See the HTML version for the code and details.

Code

subj_means_im <-
  df_i_rt |> 
  group_by(subjectid, aphantasia, congruence) |>
  reframe(rt = mean(rt)*1000)

(
  plot_models_full(subj_means_im, emmeans_im, y_min = 400, y_max = 1110, size = 4) +
  plot_models_zoomed(subj_means_im, emmeans_im, size = 4) &
  theme(
    legend.position = "top",
    legend.title = element_blank(),
    legend.text = element_text(size = 20),
    legend.margin = margin(0, 150, 25, 0),
    axis.title.y = element_text(size = 20),
    axis.text.x = element_text(size = 18)
  )
) +
  plot_layout(
    guides = "collect",
    widths = c(1, 1.2)
  )

# ggsave(here("figures/model-implicit.png"), dpi = 600)

Figure 6.1: Subject means and model-estimated means per group and condition in the implicit task.

Code

subj_means_ex <-
  df_e_rt |> 
  group_by(subjectid, aphantasia, congruence) |>
  reframe(rt = mean(rt)*1000)

(
  plot_models_full(subj_means_ex, emmeans_ex, y_min = 350, y_max = 1450, size = 4) +
  plot_models_zoomed(subj_means_ex, emmeans_ex, size = 4) &
  theme(
    legend.position = "top",
    legend.title = element_blank(),
    legend.text = element_text(size = 20),
    legend.margin = margin(0, 150, 25, 0),
    axis.title.y = element_text(size = 20),
    axis.text.x = element_text(size = 18)
  )
) +
  plot_layout(
    guides = "collect",
    widths = c(1, 1.2)
  )

# ggsave(here("figures/model-explicit.png"), dpi = 600)

Figure 6.2: Subject means and model-estimated means per group and condition in the explicit task.

Session information for reproducibility

═════════════════════════════════════════════════════════════════════════

Analyses were conducted using the R Statistical language (version 4.4.1; R Core
Team, 2024) on Windows 11 x64 (build 22631)

Packages used:

  - quarto (version 1.4.4; Allaire J, Dervieux C, 2024)
  - qqplotr (version 0.0.6; Almeida A et al., 2018)
  - lme4 (version 1.1.35.5; Bates D et al., 2015)
  - Matrix (version 1.7.0; Bates D et al., 2024)
  - effectsize (version 0.8.9; Ben-Shachar MS et al., 2020)
  - ggbeeswarm (version 0.7.2; Clarke E et al., 2023)
  - statmod (version 1.5.0; Giner G, Smyth GK, 2016)
  - simr (version 1.0.7; Green P, MacLeod CJ, 2016)
  - lubridate (version 1.9.3; Grolemund G, Wickham H, 2011)
  - ggpubr (version 0.6.0; Kassambara A, 2023)
  - labelled (version 2.13.0; Larmarange J, 2024)
  - emmeans (version 1.10.3; Lenth R, 2024)
  - parameters (version 0.22.0; Lüdecke D et al., 2020)
  - performance (version 0.12.0; Lüdecke D et al., 2021)
  - easystats (version 0.7.2; Lüdecke D et al., 2022)
  - see (version 0.8.4; Lüdecke D et al., 2021)
  - insight (version 0.20.1; Lüdecke D et al., 2019)
  - bayestestR (version 0.13.2; Makowski D et al., 2019)
  - modelbased (version 0.8.8; Makowski D et al., 2020)
  - report (version 0.5.8; Makowski D et al., 2023)
  - correlation (version 0.8.5; Makowski D et al., 2022)
  - here (version 1.0.1; Müller K, 2020)
  - tibble (version 3.2.1; Müller K, Wickham H, 2023)
  - datawizard (version 0.11.0; Patil I et al., 2022)
  - patchwork (version 1.2.0; Pedersen T, 2024)
  - R (version 4.4.1; R Core Team, 2024)
  - pacman (version 0.5.1; Rinker TW, Kurkiewicz D, 2018)
  - openxlsx (version 4.2.5.2; Schauberger P, Walker A, 2023)
  - ggplot2 (version 3.5.1; Wickham H, 2016)
  - forcats (version 1.0.0; Wickham H, 2023)
  - stringr (version 1.5.1; Wickham H, 2023)
  - tidyverse (version 2.0.0; Wickham H et al., 2019)
  - readxl (version 1.4.3; Wickham H, Bryan J, 2023)
  - dplyr (version 1.1.4; Wickham H et al., 2023)
  - purrr (version 1.0.2; Wickham H, Henry L, 2023)
  - readr (version 2.1.5; Wickham H et al., 2024)
  - scales (version 1.3.0; Wickham H et al., 2023)
  - tidyr (version 1.3.1; Wickham H et al., 2024)

═════════════════════════════════════════════════════════════════════════